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Abstract 

Most  advanced  composite  materials  have  nonlinear  stress  versus  strain  curves — grossly  so  in 
the  lamina’s  shear  direction.  If  linear  approximations  are  used,  overly  conservative  design  will 
result.  This  thesis  expands  the  usefulness  of  the  Automated  Structural  Optimization  System 
(ASTROS)  finite  element  (FE)  program  by  incorporating  a  newly  developed  capability  to  include 
material  nonlinearities.  By  dividing  a  FE  problem  into  several  intervals,  and  incorporating  the 
changing  nature  of  the  nonlinear  material  properties  (E,  and  v12)  within  each  element,  a 
good  approximation  of  the  laminate’s  actual  stresses  and  strains  can  be  obtained  without  the  undue 
experimentation  involved  in  testing  every  individual  configuration.  Each  material  property  is 
updated  in  terms  of  tension  and  compression,  as  appropriate  (G12  reacts  the  same  in  tension  and 
compression).  Updating  the  FE  model’s  material  properties  is  possible  in  the  material  axes  through 
comparison  to  cubic  spline  curves  made  from  experimental  lamina  data.  The  resulting  piecewise 
linear  analysis  accurately  captures  the  nonlinear  nature  of  the  changing  properties.  This  research 
considered  plane  stress  flat  plates  with  and  without  circular  discontinuities  (holes).  Graphite/Poly- 
etheretherketone  (GR/PEEK)  is  an  advanced  composite  material  that  has  good  high  temperature 
strength  properties  and  nonlinear  material  properties. 
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SOLVING  NONLINEAR 
COMPOSITE  MATERIAL  PROBLEMS 
USING  THE  AUTOMATED  STRUCTURAL 
OPTIMIZATION  SYSTEM  (ASTROS) 


I.  Introduction 


Motivation 

Imagine  designing  something  complex,  such  as  a  house,  without  the  different  interests  working 
together.  Fortunately,  general  contractors  coordinate  the  efforts  of  the  varying  segments  to  assure 
a  quality  product  Unknowingly,  concurrent  engineering  has  been  understood  by  the  builders  of 
houses  for  years.  It  simply  makes  sense. 

In  the  Department  of  Defense,  concurrent  engineering  describes  the  streamlining  and  overhaul¬ 
ing  of  the  design  process.  Wright  Laboratories,  at  Wright-Patterson  AFB  in  Ohio  uses  a  concurrent 
engineering  tool  for  aerospace  designers:  the  Automated  Structural  Optimization  System  (AS¬ 
TROS).  ASTROS  allows  the  smooth  and  successful  completion  of  the  roughest  of  all  design 
phases,  that  of  preliminary  design.  Figure  1  shows  that  preliminary  design  occurs  after  conceptual 
design  (33). 

Although  conceptual  design  is  usually  very  sketchy,  preliminary  design  provides  a  model  that 
will  be  within  five  percent  of  the  actual,  or  final,  design  (33).  This  means  preliminary  design 
provides  maximum  payoff  within  the  design  process  (33).  ASTROS  brings  together  mechanical, 
aeroelastic,  and  control  system  design.  Since  ASTROS  uses  finite  elements  (see  Chapter  II),  it  also 
allows  sensitivity  analysis  and  optimization.  First  delivered  in  1987,  ASTROS  is  new  and  contin¬ 
ues  to  be  improved  (18). 

Purpose 

This  present  research  affords  ASTROS  nonlinear  material  capabilities,  even  for  composites. 

For  clarity,  the  term  composites  will  be  limited  here  to  specifically  mean  more  than  just  an 
increase  of  material  properties  in  desired  directions  (which  would  include  plywood  and  reinforced 
concrete).  Here,  composite  will  denote  those  non-isotropic  materials  becoming  common  within 
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Figure  1.  Elements  of  Design 


the  aerospace  industry,  such  as  graphite  fibers  encased  in  a  matrix  of  epoxy.  Advanced  composite 
materials  (ACMs)  are  the  next  generation  of  aerospace  industry  composites — they  have  even 
further  enhanced  strength  and  stiffness  properties.  They  hold  very  promising  opportunities  for  the 
future.  Advanced  composites  either  have  metallic  or  organic  matrices.  Most  have  organic  matri¬ 
ces,  which  are  nonlinear  (20).  Graphite/Polyetheretherketone  (Gr/PEEK)  is  such  a  material.  It 
holds  special  promise  for  the  future,  because  it  adds  good  high  temperature  strength  to  the  other 
advantages  of  an  ACM.  It  may  be  used  in  projects  such  as  the  National  Aerospace  Plane,  which 
will  be  subjected  to  high  temperatures.  Gr/PEEK’s  nonlinear  properties  make  it  an  ideal  choice  for 
this  current  research.  The  term  nonlinear  means  that  the  normal  assumption  of  Hooke’s  l  aw: 

o  =  Ee 
o  =  stress 

£  =  Modulus  of  Elasticity 
£  =  strain 

is  not  valid.  This  would  be  due  to  a  non-constant  slope,  E,  of  the  stress-strain  curve.  Kinematic 
nonlinearity  is  discarded  since  most  composites  cannot  sustain  more  than  5  percent  strain  before  failure 
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(11).  More  than  being  interesting,  nonlinear  materials  hold  the  key  of  ament  and  impending  design. 
The  tailored  composite  materials  of  today  and  for  tomorrow  need  to  be  better  understood  to  be  fully 
exploited  Design  will  be  incorrect,  possibly  much  too  conservative,  when  nonlinear  materials  are 
approximated  as  linear  (20). 

Background  and  Overview 

As  nonlinear  laminae,  or  plies,  are  stacked  to  form  laminates,  errors  from  using  linear  theory 
will  compound  erratically,  and  are  not  intuitively  easily  accounted  for.  Although  there  are  many 
reports  on  the  linear  behavior  of  composite  lamina  and  their  interaction  as  laminates  (25,28,29,30) 
research  on  nonlinear  behavior  has  been  more  limited  (6,9,10,11,12,16,17,20,23,26,27). 

Several  experiments  have  been  conducted  on  the  nonlinear  properties  of  composite  materials 
(CMs)  and  ACMs.  First,  some  attempts  successfully  disclosed  the  nonlinear  laminae  properties 
from  testing  laminates,  then  others  determined  laminate  properties  from  testing  laminae.  Early 
attempts  focused  on  the  most  nonlinear  of  moduli,  the  shear  modulus,  G12. 

Petit  first  determined  a  nonlinear  lamina  G12  response  by  testing  (±45),  coupons,  assuming  no 
coupling  of  the  other  two  moduli,  Ex  and  E2,  which  were  approximated  as  linear  (20).  Hahn  and 
Tsai,  still  considering  Ex  and  E2  linear,  considered  coupling  in  their  derivation,  and  also  noted  the 
importance  of  considering  variations  in  Vi2,  Poisson’s  ratio  (11). 

Rosen  was  able  to  obtain  G12  directly  from  test  results,  improving  accuracy  by  finding  G12 
independently  of  £\  and  E2  (23). 

Lifshitz  further  refined  the  results  for  a  nonlinear  G12  by  only  requiring  axial  strains  (16). 
Transverse  strain  measurements  are  avoided  since  they  are  quite  susceptible  to  error.  G12  was  then 
calculated  from  known  values  of  Eu  E2  and  Vi2,  and  measured  values  of  axial  strain. 

Although  Hahn  and  Tsai  initially  felt  that  closed-form  representations  from  laminates  were 
impossible  (11),  work  has  progressed  in  this  area.  Hahn  was  able  to  calculate  a  nonlinear  laminate 
G12  from  lamina  data,  which  he  confirmed  experimentally  (10).  Petit  and  Waddoups  were  able  to 
produce  reasonably  good  results  for  Ex,  E2  and  Gx2  in  nonlinear  laminates  from  their  lamina 
properties  (20). 

Sandhu  has  led  the  efforts  in  analyzing  the  nonlinear  properties  of  Gr/PEEK.  Under  simple 
loading  condi tioar,  Sandhu  was  able  to  analytically  relate  laminate  behavior  to  that  in  individual 
laminae  (26).  Using  cubic  spline  interpolation  (see  Chapter  II)  within  incremental  constitutive 
equations,  he  was  able  to  determine  the  nonlinear  Eu  E2  and  GJ2  as  well  as  update  vJ2.  Daniels, 
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Fisher  and  Martin  have  confirmed  Sandhu’s  work  experimentally  with  their  respective  research  on 
the  responses  of  Gr/PEEK  to  various  loading  conditions  (6,9,17). 

This  present  research  will  incrementally  incorporate  all  three  nonlinear  moduli  while  also 
updating  Poisson’s  ratio  into  the  convenience  and  versatility  of  ASTROS.  Although  only  axial 
tension  will  be  applied,  compression  will  be  felt  in  some  elements,  as  well  in  some  material 
directions  due  to  Poisson’s  effect.  The  generated  tests  will  consider  composite  plates  loaded  in 
tension,  including  (0),  (90),  and  (±45),  laminates  using  Fisher’s  high  temperature  data. 

This  will  be  accomplished  by  the  following  procedure.  First,  ASTROS  is  run  with  nominal 
values  with  an  increment  of  the  total  load  or  displacement  to  be  applied.  The  resulting  stresses  and 
strains  for  each  element  and  angular  orientation  of  the  finite  element  model  are  compared  to 
experimental  data  for  the  moduli.  The  strain  values  update  E\ ,  £2,  G\2  and  V12  for  each  finite 
element  using  cubic  spline  interpolation  through  experimental  data  points.  After  updating  these 
properties  on  the  ASTROS  material  cards,  ASTROS  is  run  for  the  successive  increment  By 
summing  each  element’s  incremental  stresses  and  strains,  the  nonlinear  data  results. 

This  present  research  shows  the  ability  to  predict  nonlinear  laminate  results  in  different 
orientations  and  ply  lay-ups,  even  with  discontinuities,  or  holes.  Other  applications  include 
modeling  nonlinear  materials  that  are  not  composites.  The  overall  result  is  that  nonlinear  materials 
can  successfully  be  used  to  design  within  the  ASTROS  framework. 
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II.  Theory 

This  thesis  made  use  of  several  theories  and  concepts.  This  chapter’s  purpose  is  to  discuss  those 
relevant  postulates,  forming  a  basis  for  the  application  seen  in  this  research.  The  areas  of  interest 
covered  include,  sequentially:  nonlinearity,  plane  stress  elasticity,  mechanics  of  composite  mate¬ 
rials,  finite  element  mechanics,  the  ASTROS  computer  program,  the  cubic  spline,  and  the  presen¬ 
tation  of  the  sum  of  these  theories  into  the  algorithm  that  this  research  has  proposed  to 
advance — incorporating  nonlinear  material  properties  into  the  Automated  Structural  Optimization 
System  (ASTROS). 

Assume  Cartesian  coordinates  throughout  this  chapter  and  throughout  the  entire  report 

Nonlinearity 

Nonlinearities  can  be  classified  into  three  categories  (3),  including  1)  Materially  nonlinear  only; 
2)  Large  displacements,  large  rotations,  but  small  strains;  3)  Large  displacements,  large  rotations 
and  large  strains. 

This  present  research  considered  only  material  induced  nonlinearities,  or  the  first  category.  By 
using  a  sufficiently  small  increment,  the  problem  can  be  considered  piece- wise  linear,  and  more 
easily  understood  (3).  This  was  the  approach  taken  to  evaluate  the  effect  of  changing  moduli  within 
the  Graphite/PEEK  Composites  considered.  This  is  an  appropriate  techniques,  even  for  com¬ 
posites,  since  even  nonlinear  materials  will  retain  their  orthotropic  nature  (4).  This  means  that  the 
moduli  are  uncoupled,  that  is,  shear  stresses  won’t  cause  normal  straining  and  normal  stresses 
won’t  cause  shear  straining. 

Plane  Stress  Elasticity 

Volumes  could  be  and  have  been  written  about  elasticity.  For  the  puqpose  of  this  present 
research,  it  will  be  sufficient  to  encapsulate  the  theory  while  expounding  on  one  of  the  tenants  of 
two  dimensional  elasticity,  namely  the  plane-stress  assumption.  Elasticity  is  the  study  of  the 
changes  in  length  and  direction  caused  by  displacement  of  a  line  element  that  joins  two  infinitesi¬ 
mally  close  points  (24).  We  will  see  that  elasticity  is  also  behind  the  theory  of  the  finite  element 
method,  which  will  be  discussed  shortly.  Elasticity  defines  equations  that  are  basic  to  finite 
element  modeling. 

If  the  planes  initially  normal  to  an  axis  remain  normal  to  it  after  a  deformation,  including  the 
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straight  lines  initially  parallel  to  the  given  axis  remaining  parallel,  plane  deformation  is  encoun¬ 
tered  (24).  Within  plane  deformation,  plane  strain  or  plane  stress  may  result.  Plane  strain  occurs 
when  there  are  no  strains  in  the  axis  considered.  Plane  stress  involves  zero  stresses  along  the  axis 
in  question.  Consider  plane  stress  for  the  out  of  (the  horizontally)  plane  axis  X3  in  Figure  2.  Plane 
stress  would  involve: 


O33  =  *13  =  *23  =  0 

Generalized  plane  stress  involves  thin  plates  whose  stress  tensor  and  displacement  vector  are 
evaluated  in  this  horizontal  plane.  For  isotropic  (having  no  preferred  directions)  and  isothermal 
conditions,  plane  stress  can  be  represented  analytically  by: 


1  -v  0 
-v  0  0 
0  0  % 


The  plane  stress  assumption  is  well  suited  for  two  dimensional  problems.  However,  It  could  be 


invalidated  by  compression-induced  buckling  into  the  third  dimension.  Because  of  this  limitation. 


Figure  2.  Flat  Plate 
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the  present  research  is  restricted  to  axial  tension  loading  and  displacements  only.  The  plane  stress 
relations  for  a  composite  material  will  be  developed  in  the  next  section. 

Mechanics  of  Composite  Materials 

ASTROS  has  been  utilized  in  this  research  for  analysis  of  both  isotropic  and  anisotropic 
materials.  Composites  are  not  isotropic.  The  fiber  and  matrix  composites  considered  in  this 
research  are  generally  orthotropic.  As  shown  in  Figure  3,  orthotropic  means  that  along  the  fiber 
axis  there  is  symmetry.  Also  notice  in  this  figure  that  the  fiber  direction  is  aligned  with  the 
coordinate  axis.  This  is  not  required  for  a  material  to  be  orthotropic,  but  analysis  of  nonaligned 
materials  requires  a  transformation. 

There  are  two  general  approaches  to  composite  materials:  considering  micromechanical  or 
macromechanical  behavior.  Micromechanical  studies  concern  themselves  with  the  fiber  and 
matrix  materials  individually.  Macromechanical  studies  are  concerned  with  the  properties  of  the 
ply  as  a  whole.  This  is  the  interest  of  the  present  research.  Since  our  nonlinear  model  is  piece- wise 
linear,  linear  elastic  mechanics  model  the  behavior  of  our  laminates. 

Lamina  Behavior.  Within  each  linear  segment,  Hooke’s  Law  holds: 

o=£e 

e  =  Go 

In  matrix  or  tensor  notation,  Hooke’s  Law  can  be  conveniently  written: 


Figure  3.  Lamina  with  Unidirectional  Fibers  (15) 
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c,  =  CyEj  or 
Ei  =  S<jOj 

fori  &j  =1,2, 3, 4,5, 6 

where  C  is  the  stiffness  matrix  and  S  is  the  compliance  matrix. 

When  the  coordinate  axes  correspond  to  the  symmetry  axes  of  the  material,  the  matrices  can  be  written 
as: 


Ol 

C„ 

C12 

C13 

0 

0 

o' 

V 

<?2 

Cl  2 

c22 

C22 

0 

0 

0 

£2 

03 

Cl  3 

C23 

C33 

0 

0 

0 

£3 

t23 

0 

0 

0 

C44 

0 

0 

Y23 

X31 

0 

0 

0 

0 

C55 

0 

Ysi 

0 

0 

0 

0 

0 

Cm 

Yl2 

where  e„  =  engineering  strain  =  2y,; 


which  reduces  to  the  following  in  plane  stress  ( a33  =  x13  =  x23  =  0  ) 


O, 

Cn  Ci2  0 

£i 

02 

.  = 

C12  C22  0 

£2 

T12 

0  0  c« 

Yl2 

or  equivalently 


El  ’ 

5n  Sa  0 

’  1 

Oi 

e2 

■  = 

5.2  $22  0 

< 

o2 

Yl2 

0  0 

x.2 

In  the  above  equations,  C,;  and  Sy  have  been  found  to  equal: 


C„  - 


Ex 

l-v12v21 


Q2  =  C21  = 


v21£. 


1  -V,2  V21 


Cm  — ' 


1  “  vi2  V2, 

C«  -  &12 
Sn=4~ 


C  _-V12_-V2l_c 
“12  =  —  =  ~c —  =  “21 


En 


-22 
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The  stiffness  matrix,  C,  is  also  referred  to  as  Q. 

When  the  coordinate  axes  and  fiber  (or  material)  axes  are  not  aligned,  as  in  Figure  4,  a 
transformation  must  rotate  the  above  equations  into  the  coordinate  axes: 


o* 

[0i, 

012  016 

cy 

“  = 

012 

022  026 

T* 

0 16 

026  066 

Y*y 

where 


0n  =  Qn +  022*4  +  2(012  +  2066)i2oj 
022  =  Qui*  +  Q22P*  +  2(012  +  2066)l2O2 
012  =  (011  +  022  “  4066)l2O2  +  012(t4  +  O*) 

066  =  (011  +  022  ~  2012  “  2066)Po2  +  066(l'4  +  C4) 

016  “  (011  -  012  ~  2066)O31  -  (022  “  012  ~  2066)l3O 
026  ~  (011  ~  012  ~  2 064 JO3!  -  (022  “  012  ~  206fi)l'3O 

where  i  =  sinrp 
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0=cos9 


<p  =  transformation  angle 

Laminate  Behavior.  On  the  macromechanical  scale,  laminate  behavior  depends  on  the  sum¬ 
mation  of  individual  lamina  behavior.  However,  individual  plies  cannot  deform  independently,  so 
equation  (1)  must  be  satisfied  by  each  layer.  The  interaction  of  the  plies  is  accounted  for  in  the 
next  section,  within  the  finite  element’s  equivalent  stiffness  matrix.  Neighboring  plies  of  different 
orientation  usually  cause  out  of  plane  stresses,  appropriately  called  interlaminar  stresses. 

Finite  Element  Modeling  (FEM) 

The  finite  element  method  is  a  numerical  procedure  for  analyzing  structures 
and  continua.  Usually  the  problem  addressed  is  too  complicated  to  be  solved 
satisfactorily  by  classical  analytical  methods.  The  problem  may  concern 
stress  analysis,  heat  conduction,  or  any  of  several  other  areas.  The  finite 
element  procedure  produces  many  simultaneous  algebraic  equations,  which 
are  generated  and  solved  on  a  digital  computer.  Finite  element  calculations 
are  performed  on  personal  computers,  mainframes,  and  all  sizes  in  between. 

Results  are  rarely  exact.  However,  errors  are  decreased  by  processing  more 
equations,  and  results  accurate  enough  for  engineering  purposes  are  obtain¬ 
able  at  reasonable  cost  (5) 

With  that  introduction.  Cook  has  set  the  stage  for  further  discussion  of  this  technique.  The  finite 
element  method  is  used  in  designing.  FEM  can  be  used  in  the  design  of  any  conceivable 
structure — the  U.S.  Air  Force  uses  FEM  exhaustively.  Although  it  is  possible  to  use  FEM  without 
a  computer  (FEM  dates  back  to  1906),  FEM  has  come  of  age  because  of  the  vast  numerical  abilities 
of  modem  computing  equipment. 

The  main  divisions  of  FEM  include:  strain-displacement  relations,  stress-deformation  relations, 
stress-strain  relations,  stiffness  equations,  element  formulation,  and  convergence.  Although  not 
separately  addressed,  the  theory  of  elasticity  will  be  evident  in  the  following  discussion. 


Strain  Displacement  Relations.  Figure  5  shows  the  deformation  caused  by  loading  -  -  01 
and  02  become  O’  1  *  and  0’2\  Recall 


e*  =  ‘ 


Iff  T- Iff  T 
L02 


-  u-dx 
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Figure  5.  Displacement  and  Distortion  of  Differential  Lengths  dx  and  dy 


3v 

Similarly,  £,  =  — 

Engineering  shear  strain  is  defined  as  the  amount  of  change  from  a  right  angle.  Assuming  small 
angles,  P,  approximately  equals  tanp,  likewise  for  p2.  Therefore,  Cook  shows  that: 

The  above  three  equations  are  the  two  dimensional  strain  displacement  relations,  which  can  be 
written  as: 


where  U  represents  displacements. 

Stress-Deformation  Relations.  This  section  could  formally  be  entitled  “Elasticity.” 
Equilibrium,  compatibility  and  the  boundary  conditions  must  be  satisfied  within  the  finite  element 
method  for  a  correct  solution,  just  as  in  elasticity.  Equilibrium  is  satisfied  when  the  forces  from 
the  stresses  on  the  differential  element’s  edges  and  from  its  body  forces  are  in  equilibrium  (which 
is  that  the  sum  of  all  forces  equals  zero  for  a  statics  problem). 

From  Figure  6,  we  see  this  must  give  the  following: 


Our  +  +  Fm  *  0  (x  direction) 
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Figure  6.  Stresses  and  Body  Forces 


+  o>0,  +  Fy  =  0  (y  direction) 


where 


notation  denotes 


and  Fx  and  Fy  respectively  represent  the  body  forces  in  the*  and  y 


directions. 

Compatibility  is  best  explained  by  example.  If  we  consider  a  Rubik’s  cube,  we  realize  that  there 
will  not  be  gaps  or  overlaps  between  elements.  This  also  must  be  computationally  true  on  the 
infinitesimal  level. 

The  enforcement  of  boundary  conditions  translates  into  prescribed  displacements  and/or 
stresses.  Again,  boundary  conditions  are  easily  comprehended  by  example.  If  you  analyze  a 
cantilevered  beam,  the  inability  of  displacements  at  the  wall  form  a  boundary  condition  for  the 
specimen. 


Stress-Strain  Relations.  Ignoring  the  effects  of  temperature  change, 

{e}  =  [5]  {o}  or  equivalently 


W-iciW 

Where  [S]  is  called  the  compliance  tensor,  and  [C]  represents  the  stiffness  tensor,  and: 

(S]  =  [C]-» 


As  you  may  recall  from  the  previous  section,  the  plane  stress  assumption  was  written  as  a  function  of 
stress  and  strain. 
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Element  Formulation .  Stepping  back,  into  the  bigger  picture  of  finite  element  modeling,  it 
would  be  useful  to  realize  that  the  problem  is  solved  as  several  smaller  blocks  are  solved.  These 
blocks  are  interconnected,  like  the  face  of  the  Rubik’s  cube  mentioned  earlier,  except  that  there  can 
be  thousands  of  these  smaller  pieces.  They  are  interconnected  at  their  boundaries  in  terms  of  the 
constitutive  equations,  which  represent  each  block  physically. 

Before  considering  the  elements  as  a  set,  a  single  two  dimensional  element  will  be  developed. 
Although  three-sided  or  triangular  elements  were  first  used,  this  present  research  was  conducted 
using  a  quadrilateral  (rectangular)  plate  element.  This  research  could  not  be  done  with  membrane 
(i.e.  two  dimensional)  elements.  ASTROS,  as  well  as  other  finite  element  programs  require  plate 
elements  to  capture  the  thickness  which  needs  to  be  considered  for  a  composite  laminate.  The 
program  PLSTR  was  used  to  show  the  accuracy  of  the  developed  algorithm.  Although  ASTROS 
and  PLSTR  formulate  their  plate  elements  differently,  the  two  programs  gave  answers  which 
correlated  to  the  sixth  decimal  place  when  a  isotropic  material  was  considered  with  the  same 
membrane  element. 

The  plate  elements  in  ASTROS  and  PLSTR  are  formulated  with  the  same  objectives  in  mind, 
but  there  are  a  few  differences  (although  they  are  not  very  significant). 

Dr.  Sanhu’s  element  is  equivalent  to  another  type  of  ASTROS  element  They  are  formulated  by 
breaking  the  quadrilateral  into  four  triangles  after  generating  a  fictitious  node  in  the  middle  (32). 
This  type  of  element  can  then  consider  a  quadrilateral  as  four  constant  strain  triangles.  The 
polynomial  which  describes  each  triangular  element  is  of  the  form: 

U  =  ax  +  by  +  c 

It  is  complete,  meaning  that  when  the  strain  is  calculated,  the  only  terms  remaining  after  the 
partial  differentiation  are  the  constants  a  and  b,  since 


Ey 


8v 

~dy 


This  results  in  four  triangular  elements  whose  strain  remains  constant  throughout  each.  The 
fictitious  central  node  is  then  removed  to  reform  the  quadrilatural  element 
The  element  ASTROS  used  in  this  study  is  isoparametric,  with  slightly  better  stiffness  properties 
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than  the  constant  strain  elements  considered  above.  ASTROS  evaluates  the  quadrilateral  element 
using  a  displacement  equation  of  the  form: 


U  =  ax  +  by  +  cxy  +  d 

When  strains  are  found,  there  remains  a  term  which  is  either  a  function  of  x  or  y;  these  strains 
are  said  to  vary  linearly.  Results  are  found  using  numerical  integration  across  the  remaining 
variable. 

The  following  discussion  will  be  in  terms  of  this  isoparametric  four-sided  element,  although  for 
very  simple  structures  a  one  dimensional  bar  or  rod  element  may  be  used  as  shown  in  Figure  7. 
The  linear  strain  quadrilateral  plate  element  (shown  in  Figure  8)  has  eight  degrees  of  freedom 
(DOF):  each  of  the  4  nodes  can  move  in  the  x  or  y  directions.  A  displacement  in  the  x  or  horizontal 
direction  is  called  u,  whereas  a  displacement  in  the  y  direction  is  in  the  node’s  v  direction. 

Shape  functions  are  used  to  determine  the  response  of  a  DOF.  A  matrix  of  shape  functions,  [N], 
determines  how  each  displacement,  u  and  v,  varies  with  x  when  the  corresponding  DOF  has  unit 
value  and  the  other  is  zero.  IN]  describes  how  to  interpolate  each  DOF  over  the  elements. 


“il 
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Figure  7.  Structural  Elements 
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Figure  8.  QUAD  4  Element 
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Strains  can  be  found  from  the  strain-displacement  matrix,  [B],  and  the  displacement  vector: 


£x 

y*y 


=  [B] 


«i 

«2 

«3 

V;  ’ 
V2 

Vj 

V4. 


An  individual  element  stiffness  matrix,  [k],  is  now  achievable: 


[*]=J  [Bf[E\  [B]tda 


[k]  =  [B]T[E]  \B]  tA 
where: 
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1  v  0 


(in  plane  stress) 

In  composites,  [£f]  is  parallel  to  [E];  therefore  for  each  composite  element: 

[k]  is  an  equilibrium  set  of  forces  from  the  nod*s,  or  corners,  associated  with  the  activation  of 
that  DOF.  The  global  stiffness  matrix,  [K],  is  an  interlinked  combination  of  all  the  element’s  [k], 

where 

{ D}  =s  Global  DOF  matrix 
{R}  =  Global  force  matrix. 

The  discussion  of  finite  elements  would  not  be  complete  without  the  practical  application.  There 
are  several  important  concepts  worth  mentioning.  First,  th*  elements  must  be  chosen  to  have  low 
aspect  (i.e.  length  to  width)  ratios  and  the  angles  at  the  comers  of  each  element  should  not  be  too 
extreme,  but  instead  stay  around  ninety  degrees.  The  elements  must  be  small  enough  to  accurately 
model  the  actual  situation  (this  is  known  as  convergence),  yet  as  the  number  of  elements  increase, 
so  does  computation  time  and  complexity.  Each  element  and  its  nodes  are  numbered  (according 
to  the  specific  program’s  instructions)  to  keep  track  of  them,  specifically  those  representing  areas 
of  heightened  interest.  Finally,  and  most  importantly,  it  is  essential  to  be  sure  the  model  you  give 
to  the  computer  is  what  you  think  it  is.  Wrong  answers  can  be  generated  to  extreme  accuracy,  but 
they  will  still  be  wrong. 

ASTROS 

The  computer  software  base  for  this  research  is  the  Automated  Structural  Optimization  System, 
or  ASTROS.  It  is  versatile  and  powerful.  Specifically,  ASTROS’  purpose  is  to  “provide  an 
automated  tool  for  the  preliminary  design  of  aerospace  structures  which  emphasizes  the  interdis¬ 
ciplinary  requirements  of  the  design  task”  (13).  It  has  been  developed  to  meet  the  concurrent 
engineering  objectives  of  the  defense  industry:  improved,  cost-effective,  yet  quick  design. 
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Although  concurrent  engineering  may  be  considered  a  soft  science  by  some,  its  basic  concepts 
have  been  overlooked  for  so  long  that  their  application  has  been  very  beneficial.  Concurrent 
engineering  may  be  summarized  as  the  systematic  improvement  of  product  design,  manufactura¬ 
bility,  and  quality  while  reducing  time  to  market  (31).  Considering  these  advantages,  reduced  cost 
and  improved  customer  satisfaction  should  also  follow.  In  government  terms,  concurrent  engi¬ 
neering  supports  Total  Quality  Management  (TQM).  As  Turino  defines  in  Managing  Concurrent 
Engineering,  concurrent  engineering  consists  of  the  following  five  activities: 

•  Designing  for  performance 

•  Designing  for  manufacturability 

•  Designing  for  testability 

•  Designing  for  serviceability 

•  Designing  for  compliance 

As  you  can  see,  concurrent  engineering  does  more  than  include  the  design  process;  the  design 
process  is  at  the  heart  of  it  It  is  well  known  that,  “Decisions  made  during  the  design  phase  of  a 
product  have  the  biggest  impact  on  product  cost  over  the  life  of  a  product.”  (30).  Figures  9  and  10 
graphically  show  this  (31). 

The  Department  of  Defense  is  also  interested  in  fielding  equipment  faster.  Figure  1 1  shows  that 
concurrent  engineering  saves  time  by  eliminating  the  redesign  and  reverify  loops  from  the 
engineering  process. 


Impact  on  Total  Cost 

60%  TO  95%  OF  TOTAL" 
PRODUCT  COST  IS 
DETERMINED  DURING  THE 
DESIGN  PHASE  OF  A 
PRODUCT.  EARLY 
DECISIONS  HAVE  THE 
LARGEST  IMPACT. 


Total  Invested 
Cost 


TIME 


Figure  9.  Impact  of  Design  Decisions  on  Total  Product  life  Cycle  Cost 
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Figure  10.  Concept,  Design,  Review,  and  Retrofit  Costs 


Finger  and  Dixon  have  thoroughly  reviewed  mechanical  engineering  design.  They  feel  “An 
important  area  that  has  received  little  attention  to  date  is  the  creation  of  design  environments  that 
integrate  available  tools  into  a  consistent  system  to  support  the  designer.”  (7,8).  This  brief  summary 
of  concurrent  engineering  motivates  this  research  which  improves  the  capabilities  of  ASTROS,  a 
concurrent  engineering  tool  for  designers. 

ASTROS  is  designed  to  handle  aerodynamics,  optimization,  control  systems  and  the  subject  of 
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Figure  11.  Time  to  Market 
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this  research:  static  structural  analysis  (see  Figures  12  and  13).  Figure  14  shows  the  organization 
of  this  system. 

Users  interface  with  ASTROS  in  three  possible  ways.  First,  it  is  possible  to  communicate  with 
ASTROS  through  bulk  data  entries,  which  would  represent:  the  element  and  node  sequencing, 
boundary  conditions,  and  forces  and/or  displacements.  Second,  the  solution  control  commands  are 
used  to  direct  the  computer  exactly  how  to  process  the  bulk  data.  This  includes  possible  optimi¬ 
zation  strategies,  titles  for  later  data  manipulation,  control  numbers  for  boundary  condition  and 
load  cases,  and  other  options  such  as  printing  instructions.  Third,  there  is  access  into  the  program 
software,  the  Matrix  Analysis  Problem  Oriented  Language  (MAPOL). 

The  structural  analysis  portion  of  ASTROS  is  a  program  similar  to  NAS  TRAN  in  terms  of  pre- 
and  post-processor  interfaces.  Neither  have  had  a  nonlinear  materials  feature  for  composites. 

For  laminates,  ASTROS  calculates  a  material  property  matrix,  [G].  [G]  is  an  expression 
combining  all  moduli  of  each  layer,  such  that 


CONTROL  RESPONSE 


Figure  12.  ASTROS’  Engineering  Disciplines 
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This  makes  it  possible  to  consider  a  laminate  with  plies  in  several  directions  as  a  conglomerated 
ply — with  this  new  ply’s  material  axes  coinciding  with  the  stress  coordinate  system.  The  (±45), 
laminate  has  been  successfully  considered  under  this  arrangement  with  the  stress  coordinate  system 
aligned  with  the  global,  or  x  and  y,  axes.  However,  until  more  FE  programs  adopt  this  standard, 
comparing  results  outside  of  ASTROS  will  be  limited. 


Cubic  Spline 

The  cubic  spline  is  used  to  generate  continuous  function  values,  first  derivatives,  second 
derivatives  and  integrals  from  a  set  of  data  points.  This  is  done  by  mapping  the  points  to  sections 
of  cubic  curves  and  forcing  the  first  and  second  derivatives  to  be  equal  at  the  joining  of  each  of 
these  sections.  They  have  been  found  to  be  quite  superior  to  polynomials  of  interpolation  (2). 
Simply  put,  it  is  possible  to  solve  a  system  of  linked  equations,  such  as  those  generated  to  solve  a 
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Figure  14.  The  ASTROS  System  Architecture  (14) 


spline,  if  two  additional  variables  are  known.  By  setting  the  second  derivatives  of  the  first  and  last 
data  points  equal  to  a  constant,  say  zero,  the  entire  system  could  systematically  be  solved.  In  order 
to  eliminate  any  ’ffect  of  this  forced  zero  second  derivative,  an  additional  data  point  has  been 
generated.  This  additional  point  increases  the  accuracy  of  the  last  few  regions  by  allowing  a 
smooth  curve,  through  even  the  last  data  point.  No  calculations  will  be  performed  in  the  additional 
range. 

The  spline  is  used  in  this  research  to  find  the  slopes,  or  moduli,  Eh  E2  and  Gn  from  then- 
respective  stress-strain  curves. 

The  updated  Poisson’s  ratio,  vl2,  is  the  function  value,  or  zeroth  derivative  of  the  vl2  versus  et 
(longitudinal  strain)  curve.  For  £lt£2,  and  v12  separate  splines  were  used  for  tensile  and  compres¬ 
sive  conditions. 

These  curves  can  be  generated  for  any  material  from  experimental  data. 

As  a  basis  for  the  nonlinear  material  properties  updating,  a  computer  algorithm  (22)  for  a  basic 


21 


spline  was  modified  to  accomplish  this  material  properties  research  for  each  of  the  seven  possible 
curves  (see  Appendix  A).  Five  of  these  curves  plot  stress  versus  strain  and  consider  the  slope  or 
modulus.  These  five  curves  include  the  longitudinal  material  direction  in  both  tension  and 
compression,  the  transverse  material  direction  in  tension  and  compression,  and  the  shear  direction. 
The  last  two  curves  consider  the  changing  nature  of  Poisson’s  ratio  for  elements  in  tension  and 
compressioa  These  last  two  curves  plot  Poisson’s  ratio  as  a  function  of  strain  in  the  material’s 
longitudinal  direction.  The  actual  value  for  Poisson’s  ratio  is  desired  and  not  the  slope  of  these 
curves. 

A  zero  degree  composite  laminate  plate  without  a  hole  in  tension  will  only  use  three  of  these 
seven  curves  in  a  constant  stress  analysis  in  tension.  These  include  the  longtitudal  tension  curve, 
the  transverse  compression  curve  and  the  Poisson's  ratio  curve  (in  tension  as  well  since  it  is  a 
function  of  the  longitudinal  strain). 

Summary 

This  research  was  made  possible  and  is  widely  applicable  mainly  due  to  the  flexibility  of 
ASTROS’  bulk  data  deck.  By  assigning  a  material  card  to  each  element,  the  changing  nature  of  the 
moduli  and  Poisson’s  ratio  were  captured.  By  putting  in  a  composite  property,  or  PCOMP,  card 
for  each  layer  of  each  element,  a  layered  composite  was  able  to  be  considered  nonlinearly. 

An  initial  run  with  very  small  forces  or  imposed  displacements  set  the  appropriate  beginning 
material  property  values  for  the  first  increment  This  initial  iteration  was  not  included  in  the 
summation  of  resulting  stresses,  strains,  displacements  and  reaction  forces. 

This  research  considered  both  incremental  imposed  displacements  and  edge  pressures.  Forces 
on  loaded  edge  nodes  were  used  to  represent  a  pressure  or  stress  on  the  edge  of  the  plate. 
Displacements  were  imposed  directly  to  the  edge  nodes. 

Every  iteration  of  each  test  case  ran  the  spline  program  to  update  the  material  properties,  which 
in  turn  were  used  to  evaluate  the  next  increment  of  the  imposed  stress  or  displacement  The  cycle 
was  completed  as  the  new  strains  were  used  to  find  the  updated  material  property  values  from  the 
experimental  data  using  the  spline  program.  This  constituted  the  piecewise  linear  solution  to  the 
nonlinear  properties  problem. 

Recall  from  finite  elements: 

[K]  d  =  R 
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where  [*]  =  X[*L, 


Therefore: 


Incrementally: 


In  global  axes: 


[*L  =  /JJ  [Bm\[B]dzda 

vol 

where  [£>]  relates  o  to  e 

[B]  relates  e  to  the  DOF 

=  J  \  [BV  [j^  [2*W  tB]  dA 

W./  =  fJ[BlTlAijl  [B]dA 
n 

where  [Ay]  =  I  t 

ply- 1 


~Qi,  ■  F  (Eu  E2,  V12,  Gn,  sin  9,  cos  9) 
Qij~F  (£\,  E2,  Vi2,  C12) 


\K]Ad  =  AR 


Ad  =[/n_1  A/? 


dt  =  [B]  Ad 


de, 


=  [B]  Ad  (constant  thickness  due  to  plane  stress) 


In  material  coordinates: 


<fe* 

tfer 

d£y 

,  s  . 

deT 

■ = ro 

dt. 

where  [T]  = 


cos^  sin29  2  sin9  cos9 

sin^  cos29  -2  sin9  cos9 

-sinQcosQ  sin9  cos9  cos20  -  sin29 
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Exy  =  engineering  strain  =V5  yv 


Therefore: 
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V 

V 

dt  i 

£2 

=  . 

£2 

+  • 

de 2 

£l2 

total 

£12 

previous 

Vi*1n 

current  increment 

Notice  that  the  total  strain  for  each  element  in  the  material  axes  is  the  summation  of  the  current 
increment’s  strain  with  the  strains  from  the  previous  increments.  This  is  similarly  true  for  stress, 
reaction  forces  and  displacements. 

Error  results  from  the  previous  value  for  [K]  being  used  to  evaluate  the  current  strains,  stresses, 
reaction  forces  and  displacements.  To  minimize  error,  increments  must  be  kept  small. 
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HI.  Results  and  Discussion 


This  section  summarizes  the  thesis  research.  It  is  divided  into  four  logical  units  which  represent 
the  process  and  proper  procedures  for  validating  and  using  this  research’s  results. 

Finite  Element  Model  Verification 

Mesh  Size.  There  were  two  models  used  in  this  research,  although  four  were  initially  consid¬ 
ered.  Figure  IS  shows  that  one  of  these  models  does  not  contain  a  hole.  This  simple  flat  plate  was 
used  to  verify  the  FE  and  spline  programs,  as  well  as  to  easily  contrast  constant  stress  versus 
constant  displacement  loadings.  The  other  three  models  included  meshes  of  varying  coarseness, 
but  each  of  these  had  similar  overall  dimensions.  They  were  6.0  inches  long,  1 .2  inches  wide,  .084 
inches  thick  and  each  had  a  central  circular  discontinuity  of  0.4  inch  diameter.  Figures  16  through 
18  are  quarter  plates  of  these  three  models. 

To  validate  the  accuracy  of  these  models,  an  isotropic  material  was  considered  with  a  2000  PSI 
stress  applied  through  forces  on  edge  nodes  furthest  from  the  hole.  The  actual  stress  concentration 
factor  for  a  plate  of  this  geometry  is  3.46  (9).  Table  1  shows  the  stress  and  SCF  generated  by  the 
element  closest  to  the  hole  and  perpendicular  to  the  applied  force. 


Model  Number 

1 

2 

3 

Number  of  elements 

36 

66 

72 

Maximum  Stress 

10992 

6194 

6479 

SCF 

5.50 

3.10 

3.24 

Error 

59.0% 

10.4% 

6.4% 

Table  1.  Stress  Concentration  Factor  (SCF)  Comparison 


Table  1  shows  that  model  1  is  too  coarse  to  converge.  The  answers  gained  by  this  model  are  not 
valid.  Models  2  and  3  are  actually  variations  of  the  same  mesh.  They  have  been  used  to  show  that 
although  stress  concentrations  perpendicular  to  the  hole  are  best  modeled  by  the  smallest  possible 
element  there,  it  should  not  be  the  only  consideration  in  choosing  a  model.  The  mesh  is  an  important 
factor,  but  not  the  only  factor. 

Constant  Stress  Versus  Constant  Displacement  Increments.  The  rectangular  quarter 
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Figure  15.  Plate  without  Hole 
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Figure  18.  72  Element  Model 


plate  without  a  hole  is  ideally  suited  to  bring  to  light  another  important  consideration  in  model 
selection.  The  loading  conditions  must  also  be  considered  to  determine  the  accuracy  of  a  FE  model. 
Figures  19  and  20  show  that  even  in  the  model  without  a  hole  (and  therefore  without  complicated 
stress  fields)  there  are  slight  stresses  in  the  farfield  in  the  nonloaded  directions.  This  being  true  at 
the  boundary,  where  stresses  must  equal  the  boundary  conditions  and  applied  stresses,  shows  a 
limitation  of  the  constant  displacement  model.  This  makes  sense,  since  the  constant  stress  model, 
which  should  not  produce  stresses  in  the  nonloaded  axes  at  the  farfield,  has  a  different  displacement 
and  therefore  stress  situation  than  this  constant  displacement  model.  With  constant  stress,  the  edge 
region  parallel  to  the  discontinuity  moves  more  than  the  originally  parallel  edge  regions,  due  to  the 
weakness  from  the  discontinuity. 

The  constant  displacement  loading  will  produce  small  erroneous  stresses  in  the  nonloaded 
directions,  but  internal  to  the  ASTROS  model,  the  benefits  of  the  ease  of  comparison  of  the  linear 
to  the  nonlinear  case  have  made  both  constant  stress  and  constant  displacement  applications 
appropriate. 


Figure  19.  (0)  Rectangular  Plate  without  Hole  in  Transverse  Compression 
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PROGRAM  VERIFICAnON:ELEMENT  10(w/o  hole)  IN  SHEAR 


DISPLACEMENT  IN  INCHES 


Figure  20.  (0)  Rectangular  Plate  without  Hole  in  Shear 

Nonlinear  Program  Verification 

To  verify  that  the  developed  nonlinear  material  capability  worked  as  designed,  comparisons 
were  made  to  a  validated  nonlinear  material  solver.  The  first  comparison  modeled  a  16  layer  (0) 
plate  without  a  hole,  using  generated  spline  data  (see  Appendix  B).  Figures  19, 20  and  21  are  results 
of  this  comparison  for  the  nonlinear  ASTROS,  showing  the  stresses  resulting  from  the  applied 
incremental  displacements  in  the  x  direction.  Contrast  the  magnitude  of  the  nonloaded  axes  stresses 
in  Figures  19  and  20  to  that  in  Figure  21. 

To  validate  the  nonlinear  program,  this  plate  without  a  hole  and  one  with  a  hole  were  each 
compared  to  Dr.  Sandhu’s  nonlinear  PLSTR  program.  Figure  22  shows  the  comparison  between 
the  newly  developed  nonlinear  ASTROS,  linear  ASTROS  and  PLSTR.  The  linear  version  of 
ASTROS,  represented  by  the  dash-dot  line,  does  not  account  at  all  for  die  changing  modulus.  It  is 
evident  that  the  nonlinear  version  of  ASTROS  represented  by  the  solid  line  is  able  to  correctly 
model  the  nonlinear  stress-strain  curve,  which  is  modeled  by  the  dashed  line  of  the  already  proven 
program  PLSTR.  The  differences  between  nonliner  ASTROS  and  PLSTR  can  be  attributed  to  the 
differences  in  element  formulation  and  the  effect  of  increment  size  as  discussed  in  Chapter  II. 
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STRESS  IN  PSI  ^  STRESS  IN  PSI 


xlO5  PROGRAM  VER1F.:EL.  10(w/o  hole)  IN  LONG.  TENS. 


igure  21.  ( 0 )  Rectangular  Plate  without  Hole  in  Longitudinal  Tensile  Stress 


xIOM’ROGRAM  VERIPIC'ATIONHLEMRNT  10(w/o  hole)  IN  LONG.  TENS. 

(> - 1 - 1 - i - 1 - ' 

5  _  -  NONLINEAR  ASTROS 

—  - PLSTR  (DR.  SANDI1US  PROGRAM) 

,  -  LINEAR  ASTROS 


3  r 


Figure  22.  (0)  Rectangular  Plate  without  Hole  Comparison 
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A  plate  with  a  hole  comparison  completes  our  program  verification  The  following  series  of 
figures  were  ran  in  the  66  element  model  of  Figure  17.  This  model  was  used  since  ASTROS,  both 
the  linear  and  the  nonlinear  version,  had  trouble  solving  the  72  element  mesh  from  Figure  18  in  a 
composite  material  situation  (due  to  restrictions  within  ASTROS  that  have  been  separately 
addressed  with  the  ASTROS’  Program  Management  Office  personnel). 

Since  the  splines  used  both  in  Dr.  Sandhu’s  PLSTR  and  the  nonlinear  ASTROS  are  derived 
similarly,  it  is  not  surprising  to  see  in  Figures  23  through  27  that  die  slope  at  each  data  point  is  the 
same.  Clearly,  from  this  and  the  preceding  example,  the  nonlinear  version  of  ASTROS  is  working 
properly  to  incorporate  the  effects  of  changing  material  properties.  Figures  28  and  29  contrast  the 
amount  of  displacement  incurred  on  the  nonlinear  models  with  that  of  the  linear  model  at  these 
same  stress  levels.  It  is  easy  here  to  see  the  effect  of  these  varying  amounts  of  stress,  but  in  other 
examples,  a  displacement  plot  will  not  always  effectively  show  these  varying  effects  of  the  amount 
of  calculated  stress.  The  force  predicted  through  linear  approximations  that  would  be  required  to 
give  the  same  displacement  as  actually  occurred  in  the  nonlinear  material  would  be  41.1  percent 
too  high. 

The  nonlinear  ASTROS  data  is  shown  in  the  solid  lines  in  Figures  23  through  27.  It  is  evident 

PROGRAM  VERlFICATION:ELEMENT  I  IN  TRANSVERSE  COMPRESSION 


Figure  23.  Program  Verification — Element  1  in  Transverse  Compression 
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S  IN  PSI 


Figure  24.  Program  Verification — Element  6  in  Longitudinal  Tension 


Figure  25.  Program  Verification — Element  13  in  Shear 


34 


STRESS  IN  PS  I  STRESS  IN  PSI 


x  104  PROGRAM  VERIFICATION:  ELEMENT  21  IN  LONG.  TENSION 


Figure  26.  Program  Verification — Element  21  in  Lonitudinal  Tension 


Figure  27.  Program  Verification — Element  62  in  Longitudinal  Tension 
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NONLINEAR  ASTROS  16  PLY  AT  6  DEG  AT  .667  ENFORCED  DISPLACEMENT  (XI 6) 


Figure  28.  Program  Verification — Nonlinear  Displacement  Plot 
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L 1  NEAR  ASTROS  16  PLY  AT  »  DEG  LOADED  WITH  EQUIVALENT  FORCES 


Figure  29.  Program  Verification — Linear  Displacement  Plot 
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that  nonlinear  ASTROS  models  the  changing  nature  of  the  moduli.  These  figures  also  show  the 
importance  of  appropriate  increment  sizes.  The  increments  must  sufficiently  model  the  changing 
nature  of  the  moduli  or  error  in  the  amount  of  total  stress  will  accumulate.  The  error  in  total  stress 
accumulates — this  means  that  error  induced  early  is  not  corrected  out  in  the  nonlinear  ASTROS. 
The  effect  of  increment  size  will  be  discussed  in  the  following  section. 

Iteration  Size  Study 

The  results  of  the  study  of  the  effect  of  increment  size  on  accuracy  confirmed  expectations. 
Figures  30  through  34  show  that  as  increment  size  increases,  so  does  the  amount  of  error  in 
resulting  stress.  Error  starts  out  small,  but  increases  as  increment  size  is  increased.  Each  increment 
represents  a  realignment  of  the  actual  modulus  at  a  given  strain.  Notice  that  as  increments  occur  at 
the  same  value  of  displacement,  the  curves  run  parallel  to  each  other.  This  is  due  to  the  modulus, 
or  slope,  being  corrected  in  each  curve  at  that  value. 

Results  verified  several  premises.  First,  the  shape  of  the  curve  being  modeled  is  quite  important. 
Each  element  attempts  to  follow  the  spline  curve  generated  from  the  experimental  data.  Given  an 
infinite  number  of  increments,  each  element’s  curve  would  trace  the  spline  curve  exactly.  This  is 


ITERATION  STUDY  ELEMENT  I  IN  TRANSVERSE  COMP. 


Figure  30.  Effect  of  Increment  Size — Element  1 
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STRESS  IN  PSI 


x  105  ITERATION  STUDY:  ELEMENT  6  IN  LONG.  TENSION 


DISPLACEMENT  IN  INCHES  xlO3 


Figure  31.  Effect  of  Increment  Size — Element  6 


ITERATION  STUDY.  ELEMENT  13  IN  SHEAR 


Figure  32.  Effect  of  Increment  Size — Element  13 
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STRESS  IN  PS  I  STRESS  IN  PSI 


x  104  ITERATION  STUDY:  ELEMENT  2 1  IN  LONG.  TENSION 


Figure  33.  Effect  of  Increment  Size— Element  21 


X 104  ITERATION  STUDY:  ELEMENT  62  IN  LONG.  TENSION 


Figure  34.  Effect  of  Increment  Size— Element  62 
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not  practical  and  also  not  required.  Good  results  are  generated  from  approximately  twenty  well 
placed  increments.  Of  course,  the  actual  shape  of  the  curve  may  warrant  an  increased  number  of 
increments.  The  first  curve  to  depart  from  the  linear  region  must  be  catered  to.  Second,  since  error 
is  not  corrected  out,  the  rule  is  that  initial  increments  should  stay  small  to  minimize  this  error. 
Finally,  it  has  been  shown  possible  to  model  very  linear  beginnings  of  curves  in  a  few  steps — this 
is  inversely  related  to  the  premise  that  a  quickly  changing  slope  requires  more  increments  to  keep 
stress  values  accurate. 

Linear  Venus  Nonlinear  Analysis 

The  finite  element  model  has  been  found  to  be  acceptable,  and  this  nonlinear  version  of 
ASTROS  has  been  shown  to  accurately  model  a  nonlinear  material.  This  section  considers  actual 
Gr/PEEK  experimental  data  for  the  spline.  These  curves  are  shown  in  Appendix  C.  This  section 
shows  the  results  of  linear  ASTROS  compared  to  the  newly  developed  nonlinear  ASTROS  for  (0), 
(±45)„  and  (90)  laminates  with  both  constant  displacement  and  constant  stress.  Gr/PEEK  is 
nonlinear  mainly  in  the  transverse  and  shear  directions. 

The  ( 0 )  Laminate.  This  laminate  shows  the  most  apparent  nonlinearities.  This  makes  sense 
since  the  laminate  was  loaded  in  the  x  direction,  corresponding  to  the  material’s  longitudinal 
direction — which  had  the  most  strength.  This  allowed  for  the  nonlinearities  to  develop.  Figures  35 
through  49  show  the  developing  nonlinear  characteristics  of  the  Gr/PEEK.  The  first  set.  Figures  35 
from  39,  has  a  constant  displacement  (in  the  x  direction),  whereas  Figures  40  through  44  are  in 
constant  applied  stress.  Figures  45  through  49  show  that  these  can  be  compared  by  considering  the 
average  displacement  of  the  constant  stress  model.  Figure  35  shows  the  effects  of  the  decreasing 
transverse  compressive  modulus.  Figures  36,  38  and  39  show  the  slight  nonlinearties  of  the 
longitudinal  modulus  in  tension.  The  shear  of  Figure  37  shows  nonlinearities  beginning  to  develop. 
Notice  that  not  all  nonlinearities  will  result  in  a  decrease  of  the  linearly  predicted  stress.  Since  the 
slightly  nonlinear  Gr/PEEK  longitudinal  stress  versus  longitudinal  strain  in  tension  curve  actually 
increases  slightly  in  slope,  the  actual  or  nonlinear  stress  results  are  greater  than  the  linear  in  constant 
displacement  This  reveals  that  the  modulus  is  able  to  increase  as  well  as  decrease.  This  occurs  as 
the  early  section  of  the  spline  curve  arches  steeply.  The  spline  curve’s  modulus  starts  equivalent 
to  the  linear  modulus,  but  increases  before  decreasing.  The  above  comments  also  apply  to  the 
corresponding  constant  applied  stress  model’s  figures  (Figures  40  through  44).  The  correlation  is 
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STRESS  IN  PSI  STRESS  IN  PSI 


Figure  36.  (0)  Element  6  in  Longitudinal  Tension 
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STRESS  IN  PS  I 


Figure  39.  (0)  Element  65  in  Longitudinal  Tension 


Figure  40.  (0)  Element  l  in  Transverse  Compression — Constant  Applied  Stress 
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STRESS  IN  PSI 


Figure  41.  (0)  Element  6  in  Longitudinal  Tension — Constant  Applied  Stress 


Figure  42.  (0)  Element  13  in  Shear— Constant  Applied  Stress 
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STRESS  IN  ESI 


Figure  44.  (0)  Element  65  in  Longitudinal  Tension — Constant  Applied  Stress 


46 


Figure  45.  (0)  Element  1  in  Transverse  Compression — Comparison 


Figure  46.  (0)  Element  6  in  Longitudinal  Tension — Comparison 
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STRESS  IN  PSI  |  STRESS  IN  PSI 


xlO4  (0)  ELEMENT  38  IN  LONG.  TENSION  -  COMPARISON 


Figure  48.  (0)  Element  38  in  Longitudinal  Tension — Comparison 
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X I0-*  (0)  ELEMENT  65  IN  LONG.  TENSION  -  COMPARISON 


Figure  49.  (0)  Element  65  in  Longitudinal  Tension — Comparison 


shown  in  Figures  45  through  49.  Notice  the  error  that  would  be  expected  from  linear  approxima¬ 
tions.  This  highlights  the  shortcomings  of  linear  analysis  and  the  need  for  nonlinear  analysis  of 
nonlinear  materials.  Figures  50  and  5 1  show  the  differences  between  displacements  of  the  constant 
applied  stress  and  the  constant  displacement  models. 

The  ( 90)  Laminate.  The  (90)  laminate  breaks  quickly  in  tensile  loading  since  only  the  weak 
matrix  (PEEK)  withstands  the  applied  displacements.  Figures  52  through  56  show  a  comparison 
to  those  elements  discussed  in  the  previous  cases.  Because  this  (90)  laminate  fails  quickly,  or  with 
much  less  developed  stress,  the  nonlinear  curves  are  closer  to  the  linear  approximations.  The  (90) 
curves  have  been  presented  to  enable  comparison  with  the  (0)  case.  For  example,  this  is  accom¬ 
plished  by  plotting  the  stresses  of  Figure  17’s  element  6  longitudinally  for  the  (0)  laminate  and 
transversely  for  the  (90)  laminate.  Therefore,  both  are  along  die  global  x  axis.  However,  the 
longitudinal  and  transverse  stress-strain  curves  are  different,  so  that  the  comparisons  are  different. 

The  (±45)$  Laminate.  None  of  the  moduli  for  the  (±45)*  laminate  showed  excessive  nonlin¬ 
earity  before  failure  began.  Figures  57  through  61  show  the  analysis  done  with  increments  of 
displacement 
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STRESS  IN  PSI  STRESS  IN  PSI 


Figure  52.  (90)  Element  1  in  Longitudinal  Compression 


(90)  ELEMENT  6  IN  TRANSVERSE  TENSION 


Figure  53.  (90)  Element  6  in  Transverse  Tension 
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RESS  IN  PSI  STRESS  IN  PSI 


Figure  54.  (90)  Element  13  in  Shear 


Figure  55.  (90)  Element  38  in  Transverse  Tension 
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STRESS  IN  PSI 


Figure  56.  (90)  Element  65  in  Transverse  Tension 


Figure  57.  (±45)  Element  1  in  Transverse  Compression 
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STRESS  IN  PS  I 


Figure  59.  (±45)  Element  13  in  Shear 
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NPSI 
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Recall  that  ASTROS  models  a  layered  composite  as  a  summation  of  each  ply’s  material 
properties.  Figures  57  through  61  show  that  nonlinear  materials  can  be  successfully  incorporated 
within  ASTROS.  Recall  that  the  material  axes  for  this  combined  laminate  is  along  the  global  (the 
x  and  y  directions)  axes. 

Note  in  each  of  Figures  57  through  61  that  for  a  (±45),  laminate,  the  combination  of  the  445 
fiber’s  materia]  properties  with  the  -45  fiber’s  material  properties  produce  a  laminate  whose 
nonlinear  effects  cancel  each  other  out. 

Conclusion 

The  effect  of  increment  size  on  this  nonlinear  algorithm  has  revealed  that  for  reasonable 
accuracy,  at  least  20  increments  should  be  considered.  This  of  course  depends  on  the  second 
derivation  of  the  spline’s  stress-strain  curves.  More  linear  sections  of  these  curves  may  be 
evaluated  with  less  increments,  whereas  severely  nonlinear  curves  must  have  many  increments. 
The  most  nonlinear  curve  analyzed  must  be  catered  to,  demanding  several  small  increments. 

By  analyzing  both  constant  displacements  and  constant  applied  stresses,  this  nonlinear  materials 
capability  allows  ASTROS  users  flexibility. 

There  is  benefit  in  modeling  materials  nonlinearly;  as  materials  get  even  more  nonlinear  than 
they  currently  are,  this  benefit  increases.  Further  research  could  consider  adding  failure  conditions 
to  this  nonlinear  finite  element  model.  This  research  has  accomplished  what  it  has  proposed  to  do. 
It  has  been  possible  to  model  nonlinear  materials,  including  composite  laminates,  on  the  automated 
structural  optimization  system  (ASTROS)  and  expand  the  capabilities  of  this  FE  analysis  tool. 
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Appendix  A 


cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

SPLINE  FOR  THESIS 

THIS  PROGRAM  IS  PART  OP  CAPTAIN  BOB  BOLHA’S  THESIS. 

IT  RECEIVES  STRAIN,  STRESS,  REACTION  FORCE,  AND 
DISPLACEMENT  INFORMATION  FROM  ASTBOS,  A  FINITE  ELEMENT 
PROGRAM.  THIS  INPUT  DATA  IS  COMBINED  WITH  EXPERIMENTAL 
DATA  ON  MATERIAL  PROPERTIES,  TO  DETERMINE  WHEAE  ON  THE 
STRESS  Vs  STRAIN  CURVE  THE  ELEMENT  IS  AND  WHAT  THE 
CURRENT  VALUE  OF  THE  MODULUS,  OR  SLOPE,  IS.  THESE  MODULI 
ALONG  WITH  UPDATED  POISSON’S  RATIO  INFORMATION  IS  THEN 
OUTPUT,  ALONG  WITH  THE  SUM  OF  THE  STRAIN,  STRESS,  REACTION 
FORCES  AND  DISPLACEMENTS.  THIS  OUTPUT  IS  FED  INTO  ASTROS 
TO  EVALUATE  THE  NEXT  LOAD  INCREMENT.  THE  FINAL  OUTPUT 
ACCOUNTS  FOR  THE  NONLINEAR  MATERIAL  PROPERTIES. 

PROGRAM  SPLINE 

DIMENSION  B ( 10 , 500) , C ( 10 , 500) , DELSQY (10 , 500) , DBLY (10 , 500) , 
ft  SS2(10,500) ,S2(10,500) ,S3(10,500) ,T(10,500) ,F0RCE(5) , 

ft  Y(10,500) ,F0RCET0T(5) ,IN(10) , INI (10) ,W(10) ,HT1 (10) , 

ft  HT2(10) ,PROD(10) , DELSQS(IO) ,TNOW(10,500) ,TTOT(10,500) , 

ft  X(10,500) , H(10,S00) ,H2(10,500) ,SS(10,500) ,SS1(10,500) , 

ft  STRESNOW ( 10 , 500) ,STRESTCT(1 0,500) .DISPX(IOO) ,DISPY(100) , 

ft  DISPZ(IOO) ,DISPTOTX(100) ,DISPTOTY(100) ,DISPTOTZ(iOO) 

OPEN (UNIT= 1 1 , FILE=  ’ SPL1 5E1T . DAT ’ , STATUS= ’ OLD ’ ) 

OPEN (UNIT=12 , FILB=  ’ SPL15E2T . DAT ’ , STATUS= ’ OLD ’ ) 

OPEN (UNIT=13 , FILE= ' SPL15G . DAT ’ , STATUS= ’ OLD ’ ) 

OPEN (UNTT=14, FILE=  ’ SPL15NUT . DAT ’ , STATUS=  ’ OLD ’ ) 

OPEN (UNIT=1 5 , FILE= ’SPL15E1C.DAT’, STATUS* ’ OLD 1 ) 

OPEN (UNIT= 1 6 , FILE*  ’ SPL1 5E2C . DAT ’ , STATUS* ’ OLD ’ ) 

OPEN (UNIT=17 , FILE*  ’ SPL15GX2 . DAT ’ , STATUS* ’ OLD ’ ) 

OPEN (UNIT* 1 8 , FILE= ’ SPL 1 5NUC . DAT ’ , STATUS* ’ OLD ’ ) 

OPEN (UNIT=19 , FILE* ’SPLINE170 . DAT ’ , STATUS= ’ OLD ’ ) 

OPEN (UNIT=20 , FILE* ’ INPUTCHK . OPT ’ , STATUS= ’ OLD ’ ) 

OPEN (UNIT=21 , FILE* ’ SPL15E1T . OPT ’ , STATUS* ’ OLD 1 ) 

OPEN (UNIT=22 , FILE=  ’ SPL15E2T . OPT ’ , STATUS= ’ 01 D ’ ) 

OPEN (UNIT=23 , FILE*  ’ SPL1 5G . OPT 1 , STATUS= ’ OLD ’ ) 

OPEN (UNIT- 24 , FILE=  ’ SPL15NUT . OPT ’ , STATUS= ’ OLD ’ ) 

OPEN (UNIT=25 , FILE= ’ SPL 1 5E1 C . OPT ’ , STATUS= ’ OLD ’ ) 

OPEN (UNIT=26 , FILB=  ’ SPL15E2C . OPT ’ , STATUS= ’ OLD ’ ) 

OPEN (UNIT=27 , FILE=  ’ SPL1 5GX2 . OPT ’ , STATUS= ’ OLD ’ ) 

OPEN (UNIT=28 , FILB=  * SPL15NUC . OPT ’ , STATUS= ’ OLD  ’ ) 

OPEN (UNIT=29 ,FILE= * SPLINE170 . OPT ’ , STATUS= ’ OLD ’ ) 

OPEN (UNIT=30 ,FILE=  ’ SPT0T170 . OPT ’ , STATUS= ’ OLD ’ ) 

OPEN (UNIT*31 , FILB=  ’ SPSTRSTOT . 10 ’ , STATUS= ’ OLD  ’ ) 

EPSLN=0.01 
N0DES=89 

C  M  IS  THE  NUMBER  OF  MAT8  CARDS  GENERATED 
M=66 

C  RHO  IS  MASS  DENSITY,  IT  IS  NEEDED  FOR  MATS  CARDS 
RHO=0.056 

WRITE(6, «) ’  IF  THIS  IS  THE  1ST  ITERATION  YOU  WILL’ 

WRITE (6,*)’  NEED  TO  RUN  FRONT  HALF  TO  MAKE’ 

WRITE(6, •) ’  THE  TOTAL  STRESS  ft  STRAINS  ZERO  ->  TYPE  1  (ELSE  2)’ 
READ(5, «)  ISTART 
IF (1ST ART  NE.l)  GOTO  87 
IC0UNT=O 
DO  101  K9  =  1,4 
101  FORCETOT (K9) =0 . 0 
WRITE(31 , *)  ICOUNT 
DO  95  KS=1,M 
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DO  93  K4=l ,3 

STRESTOT (K4 ,X5) =0 . 0 
TT0T(K4,K5)=O.O 
93  CONTINUE 

WRITE (3 1 , • ) K5 , STRESTOT ( 1 , K5) , STRESTOT (2 , K5) , STRESTOT (3 , KS) 
WRITE(31 , •)  TT0T(1,K5) ,TT0T(2,K5) ,TT0T(3,K5) 

95  CONTINUE 

WRITE(31 , .)F0RCET0T(1) ,P0RCET0T(2) ,F0RCET0T(3) ,F0RCET0T(4) 

DO  97  K0=1 ,  NODES 
DISPT0TX(K0)=O.O 
DISPTOTY (KO) =0 . 0 
DISPTOTZ (KO) =0 . 0 

WRITE(31,199O)DISPT0TX(KO) .DISPTOTY (KO) ,DISPTOT2(KO) 

97  CONTINUE 

CLOSE (UNIT=31) 

OPEN (UNIT=31 , FILE= ’ SPSTRSTOT . 10 ’ , STATUS= ’ OLD ’ ) 

87  CONTINUE 
1990  F0RUAT(3E15.5) 

C 

C  THIS  PART  OF  THE  PROCRAX  LOADS  IN  THE  SPLINE  CURVE  DATA: 

DO  90  IW  =  1,8 

C  IN  IS  AN  ARRAY  TELLING  HOW  MANY  DATA  PTS  EACH  SPLINE  CURVE  HAS 
READ (10  ♦  IW,.)  IN (IK) 

WRITE (20, *)  ’ 

WRITE(20,  •)  'FOR  INPUT  FILE  \IW  .10,’:  ((SPLINE  PTS=’,IN(IW) 
WRITE (20, •) '  X  Y’ 

WRITE(6 , «)  ’IN(’ ,IW, ’)=' ,IN(IW) 

DO  70  1=1 ,IN(IW) 

READ(10  ♦  IW,.)  X ( IW , I )  ,Y(IW,I) 

WRITE(20, .)  X(IW,I),Y(IW,I) 

70  CONTINUE 

X(IW,IN(IW).1)=X(IW,IN(IW)).(X(IW,IN(IW))-X(IW,IN(IW)-1)) 

Y(IW,IN(IW)«1)=Y(IW,IN(IW))*(Y(IW,IN(IW))-Y(IW,IN(IW)-1)) 

IN(IW)=IN(IW).l 

2  INI (IW) =IN(IW) -1 

3  DO  51  1=1, INI (IW) 

H(IW,I)=X(IW,I-*1)  -X(IW,I) 
DELY(IW,I)=(Y(IW,I.1)-Y(IW,I))/H(IW,I) 

51  CONTINUE 

4  DO  52  1=2, INI (IW) 

H2(IW,I)=H(IW,I-1)«H(IW,I) 

B(IW,I)=.5«H(IW,I-1)/H2(IW,I) 

DELSQY (IW , I) = (DELY (IW, I) -DELY (IW, I- 1 ) ) /H2 (IW, I) 

S2 (IW , I) =2 . »DELSQY (TW , I) 

C (IW , I) =3 . .DELSQY (IW, I) 

52  CONTINUE 
S2(IW, 1)=0. 

S2(IW,  IN(IW)  )=0 . 

0UEGA=1. 0717968 

5  ETA=0 . 

6  DO  10  1=2, INI (IW) 

7  W(IW)=(C(IW,I)-B(IW,I).S2(IW,I-l)-(.5-B(IW,I)).S2(IW,I.l)- 
k  S2(IW,I)).0MEGA 

8  IF(ABS(W(IW))-ETA)  10,10,9 

9  ETA=ABS(W(IW)) 

10  S2(IW,I)=S2(IW,I).W(IW) 

13  IF (ETA-EPSLN) 14,5,5 

14  DO  53  1=1, INI (IW) 

53  S3(IW,I)=(S2(IW,I*1)-S2(IW,I))/B(IW,I) 

90  CONTINUE 
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WRITE(20, •)  ’ ICOUNT , K , K3 , J ,  STRESTOT (K3 , J) , STRESNOW (K3 , J) ’ , 
l  \  TOTSTRAIN  (K3,J)’ 

WRITE(30,*)’  J,  SS1(1,J),  SSI  (2,  J) ,  SS1(3,  J)  ,  SS(4,J)* 

*RITE(30,»)  ’TTOT(l,J) ,  TT0T(2,J),  TTOT(3,J)’ 

WRITE(30, •) ’TNOW(l , J) ,  TN0W(2,J),  TN0W(3,J)’ 

READ(31 , »)  ICOUNT 
C 

C  THIS  PART  OF  THE  PROGRAM  UPDATES  THE  MATERIAL  PROPERTIES: 

15  DO  61  J=1,M 

READ(10,*)  TNOW(l ,  J) ,  TN0W(2,J),  TN0W(3,J) 

READ(19, *)  STRESNOW(l.J) ,  STR£SN0W(2, J) ,  STRESN0W(3, J) 

TNOW (4 , J) =TNOW ( 1 , J) 

STRESNOW (4 , J) =STRESNOW ( 1 , J) 

READ(31 , »)  JUNK, STRESTOT (1,J).  STRESTOT (2 , J) , STRESTOT (3 , J) 

READ (31  ,  ♦)  TTOT(l,J) ,  TTOT(2,J),  TT0T(3,J) 

STRESTOT (4 , J) =STRESTOT ( 1 , J) 

TTOT (4 , J) =TTOT(l , J) 

DO  89  K3=l,4 

STRESTOT (K3,J)=STREST0T(K3,J)  ♦  STRESNOW (K3,J) 

TTOT (K3 , J)  =TTOT (K3 , J)  -  TN0W(K3,J) 

T  (K3 , J) =ABS (TTOT (K3 , J) ) 

IF  (STRESTOT (K3,J).LT.  0.  )  k=k3*4 
IF  ( STRESTOT (K3,J).EQ.  0.  )  k=k3 
IF  (STRESTOT (K3,J).GT.  0.  )  k=k3 

WRITE (20 ,  * )  ICOUNT , K , K3 , J , STRESTOT (K3 , J) , STRESNOW (K3 , J) , TTOT (K3 , J) 

16  1=1 

54  IF(T(K3,J)-X(K,1))  58,17,55 

55  IF  (T  (K3 ,  J)  -X  (K ,  IN  (K)  - 1 ) )  57,59,58 

56  IF (T (K3 , J) -X (K , I) )  60,17,57 

57  1=1-1 
GOTO  56 

58  PRINT  44, J 

44  FORMAT (14, 24HTH  ARGUMENT  OUT  OF  RANGE) 

WRITE(6,.)  J.’ARG  OUT  OF  RANGE  IN  COL’, K, ’ (IF  C0L>4,  C)=’,T(K3,J) 
WRITE (20-K, »)  J, ’ARGUMENT  OUT  OF  RANGE  IN  COL’, K,’ (IF  C0L>4,  C) ’ 
WRITE(29 , *)  J.’ARG  OUT  OF  RANGE  IN  COL’, K, ’(IF  C0L>4,  C)=’,T(K3,J) 
WRITE(30, *)  J, 'ARGUMENT  OUT  OF  RANGE  IN  COL’, K, ’(IF  C0L>4,  C)  ’ 

GOTO  89 

59  I=IN(K)  -1 

60  1=1-1 

17  HT1 (K) =T (K3 , J) -X (K , I) 

HT2  (k)  =t  (k3 ,  J)  -X  (K ,  1+ 1) 

PROD (K) =HT1 (K) *HT2 (K) 

SS2  (K ,  J)  =S2  (K ,  I)  -HT1  (K)  «S3  (K ,  I) 

DELSQS (K) =  (S2 (X , I) +S2 (K , 1+ 1) *SS2 (K , J) ) /6 . 

SS  (K3 ,  J)  =Y  (K ,  I)  -HT1  (K)  .DELY  (K ,  I)  +PROD  (K)  .DELSIJS  (K) 

SSI (K3 , J) =DELY (K , I) ♦ (HT1 (K) *HT2 (K) ) *DELSQS (K)  *PROD (K) *S3 (K, I) /6 . 
IF(K  .EQ.  4)  WRITE (24,*)  SS(K, J) ,K,J, ICOUNT, STRESTOT (K, J) , 
k  STRESNOW (K , J) ,  ’  <-SS (K , J) , K , J , ICOUNT , STRESTOT (K , J) , STRESNOW (K , J) ’ 
IF (K . EQ .  8)  WRITE(28, *)  SS(X, J) ,K,J, ICOUNT, STRESTOT (K,J) , 
t  STRESNOW (K , J) , ’ < -SS (K , J) , K , J , ICOUNT , STRESTOT (K , J) , STRESNOW (K ,  J)  ’ 
IF(K.LT.4)  WRITE(20-K, *)  SSI (K, J) ,K,J, ICOUNT, STRESTOT (K,J) , 
k  STRESNOW (K,J), 'SSI, K,J,ICT, STRESTOT, NOW’ 

EF((K  .LT.  8) . AND . (K  .GT.  4))  WRITE(20*K, *)  SSI (K,J) ,K,J, ICOUNT, 
t  STRESTOT (K , J) , STRESNOW (K,J),’SS1,K,J,ICT, STRESTOT , NOW ’ 

WRITE (6 , * )  ’ STRESTOTiNOW ( ’ , k3 , ’ , ’ , j , ’ ) = ’ , STRESTOT (k3 , j ) , 
t  STRESNOW (K3.J) 

89  CONTINUE 

WRITE (29, 2000) J,SS1(1,J) ,SS1(2,J) ,SS(4, J) ,J 
2000  FORMAT (5HMAT8* , 17X, 12 , 3E16 .7 ,3H*C0 , 12) 
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WRITE (29 , 2010) J , SSI (3 , J) , RHO 
2010  FORilAT  (3H»C0 , 12 , 3X , El 6 . 7 , 32X ,  El  6 . 7) 

WRITE(30,») ’  1 ,J,SS1 (1, J) , SSI (2, J) ,SS1(3,J) ,SS(4,J) 

WRITE (30 , . ) TTOT ( 1 , J) , TTOT (2 , J) , TTOT (3 , J) 

WRITE (30,.)  TN0W(1,J) ,TN0W(2,J) ,TN0W(3, J) 

61  CONTINUE 

READ (31 , • )  FORCETOT (1 ) , FORCETOT (2) , FORCETOT ( 3) , FORCETOT (4) 
READ(19,.)  FORCE (1)  ,F0RCE(2)  ,F0RCE(3)  ,F0RCE(4) 

FORCETOT ( 1 ) =FORCETOT ( 1 )  -FORCE  ( 1 ) 

FORCETOT (2) =FORCETOT (2)  -FORCE (2) 

FORCETOT (3) =FORCETOT (3)  -FORCE (3) 

FORCETOT  (4) =FORCETOT (4)  -FORCE (4) 

DO  99  K6  =1 .NODES 

READ (31 , 1990) DISPTOTX (K6) , DISPTOTY (K6) , DISPTOTZ (K6) 

READ (19, 1990) DISPX(K6) ,DISPY(K6) ,DISPZ(K6) 

DISPTOTX (K6)=DISPT0TX(K6)  *  DISPX(K6) 

DISPTOTY (K6)=DISPT0TY(K6)  -  DISPY  (K6) 

DISPTOTZ (K6)=DISPT0TZ(K8)  -  DISPZ(K6) 

99  CONTINUE 

CLOSE (UNIT=31) 

OPEN (UNIT=31 , FILE= ’ SPSTRSTOT . 10 ’ , STATUS= ’ OLD ’ ) 

ICOUNT=  ICOUNT  -  1 
WRITE(31 , *)  ICOUNT 
DO  75  K7=l , U 

WRITE (31 ,.) K7 , STRESTOT ( 1 , K7 ) , STRESTOT (2 ,K7) , STRESTOT (3 , K7) 
WRITE (31 , .)  TT0T(1 ,K7) ,TT0T(2,K7) ,TT0T(3,K7) 

75  CONTINUE 

WRITE(31 , •)  FORCETOT (1) ,F0RCET0T(2) .FORCETOT (3) ,F0RCET0T(4) 

DO  98  Kl=l, NODES 

98  WRITE(31 , 1990) DISPTOTX (Kl) , DISPTOTY (Kl) .DISPTOTZ(Kl) 

62  CONTINUE 
CLOSE  (UNIT=U) 

CLOSE (UNIT=12) 

CLOSE (UNIT=13) 

CLOSE (UNIT=14) 

CLOSE (UNIT=15) 

CLOSE (UNIT=16) 

CLOSE (UNIT=17) 

CLOSE (UNIT=18) 

CLOSE (UNIT=19) 

STOP 

END 
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LONGITUDINAL  STRESS  (PS I) 


Appendix  B 


x  10s  LONGITUDINAL  SPLINE  CURVE  IN  TENSION 


Figure  B-l  01  (ksi)  versus  ei  ("/*.)  (tension) 
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LONGITUDINAL  STRESS  (PSD 


xlO5  LONGITUDINAL  SPLINE  CURVE  IN  COMPRESSION 


Figure  B-2  ci  (ksi)  versus  ei("4i)  (compression) 
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TRANSVERSE  STRESS  (PSI) 


R  STRESS  (PSI) 


SHEAR  SPLINE  CURVE 


Figure  B-5  T12  (ksi)  versus  yi2(±45) 


LONGITUDINAL  STRESS  (PSI) 


Appendix  C 


x  105  GR/PEEK  LONGITUDINAL  SPLINE  CURVE  IN  TENSION 


Figure  C-l  Gr/PEEK  oi  (ksi)  versus  ei(*54>)  ( tension ) 
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xlO4  GR/PEEK  LONGITUDINAL  SPLINE  CURVE  IN  COMPRESSION 


Figure  C-2  Gr/PEEK  01  (ksi)  versus  Ei("/*.)  ( compression ) 
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NSVERSE  STRESS  (PSI) 


TRANSVERSE  STRESS  (PSI) 


xlO4  GR/PEEK  TRANSVERSE  SPLINE  CURVE  IN  COMPRESSION 


Figure  C-4  Gr/PEEK  02  (ksi)  versus  £2 (*/*)  (compression) 
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SHEAR  STRESS  (PSI) 


xlO4  GR/PEEK  SHEAR  SPLINE  CURVE 


Figure  C-5  GrfPEEK  T12  (ksi)  versus  Yl2(±45) 
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POISSONS  RATIO  (NU) 


0  0.002  0.004  0.006  0.008  0.01  0.012 

LONGITUDINAL  STRAIN  (IN/IN) 


Figure  C-6  Gr/PEEK \n  (ksi)  versus  £i  (tension) 


POISSONS  RATIO  (Nil) 


LONGITUDINAL  STRAIN  (IN/IN)  xlO3 


Figure  C-7  Gr/PEEK  V12  (ksi)  versus  ei  (compression) 
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